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Abstract 

Superbursts are rare and energetic thermonuclear carbon flashes observed to occur on accreting 
neutron stars. We create the first multi-zone models of series of superbursts using a stellar evolution 
code. We self-consistently build up the fuel layer at different rates, spanning the entire range of 
observed mass accretion rates for superbursters. For all models light curves are presented. They 
generally exhibit a shock breakout, a precursor burst due to shock heating, and a two-component 
power-law decay. Shock heating alone is sufficient for a bright precursor, that follows the shock 
breakout on a short dynamical time scale due to the fall-back of expanded layers. Models at the 
highest accretion rates, however, lack a shock breakout, precursor, and the first power law decay 
component. The ashes of the superburst that form the outer crust are predominantly composed of 
iron, but a superburst leaves a silicon-rich layer behind in which the next one ignites. Comparing the 
model light curves to an observed superburst from 4U 1636-53, we find for our accretion composition 
the best agreement with a model at three times the observed accretion rate. We study the dependence 
on crustal heating of observables such as the recurrence time and the decay time scale. It remains 
difficult, however, to constrain crustal heating, if there is no good match with the observed accretion 
rate, as we see for 4U 1636-53. 

Subject headings: accretion, accretion disks — methods: numerical — nuclear reactions, nucleosyn- 
thesis, abundances — stars: neutron — X-rays: binaries — X-rays: bursts 



1. INTRODUCTION 

X-ray flares have been observed from accreting neu- 
tron stars that are similar to Type I X-ray bursts, 
but that are a thousand times more e nergetic and 
last up to a day. Normal burst s (e.g., iLewin et alJ 
119931 : IStrohmaver fe Bildstenl l2006f ) result from hydro- 
gen and helium burning to carbon and, through a 
series of alpha captures, the ap-process, and pro- 
ton captures , the rp-pro cess, to heavier elements 
(|Schatz et al.l 120011 . l2003bD . The long flares, named 
'superbursts', are attributed to the runaway thermonu- 
clear burning of carbon in a 100 m thick layer of 
ashes of normal bursts (jCumming fc Bildsten 120011 : 
IStrohmaver fc Brownl 120021 ). The day long decay is ex- 
plained as bei ng the cooling time scale o f a layer of 
that thickness (|Cumming fc Macbeth! 120041) . Because it 
takes typically about one year to build up this 100 me- 
ter thick layer, superbursts are much rarer than regular 
bursts. The first superbursts have been discovered rela- 
tively recently: it was only after the launch of the Bep- 
poSAX and RXTE observatories that enough exposure 
time was collected to be able to detect such rare e vents 
(jCornelisse et al.ll2000t IStrohmaver fc Brownl 120021 ). At 
the time of writing, 20 (candi date) superbursts have be en 
observed from 11 s ources (seelKeek fc in 't Zan d 2008 for 
an overview, a nd iKuulkerd 120091 iChenevez et al.ll2011l 
lin't Zand et al.ll2011l for recent discoveries). 

Two types of superbursts are discerned based on 
the composition of the material that is accreted from 
the companion star. Most superbursters are thought 
to accrete hydrogen-rich material. Their superbursts 
are energetic, but the peak brightness does not reach 
the Eddington limit. S uperbursts have been ob- 
served from 4U 0614+91 (jKuulkers et alJ l2010f) and 
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4U 1820-30 (IStrohmaver fc Brownll2002l and a candidate 
lin't Zand et al.H2011l ). These sources are so-called ultra- 
compact X-ray binaries (UCXBs). UCXBs have a bi- 
nary period of less than 80 minutes. In such a small 
orbit, stable mass transfer by Roche-lobe overflow can 
only occur from an evolved star that has lost its hydro- 
gen envelope. The material accreted onto the neutron 
star, therefore, contains no hydrogen, but may contain 
helium. The superburst from 4U 1820 — 30 reached the 
Eddington luminosity and displayed photospheric radius 
expansion. For 4U 0614+91 the onset of the superburst 
was not observed. 

Most superbursting sources have a high accretion 
rate M of at least 10% of the Eddington-limited rate 
-^Edd at the time of the superburst. Exceptions are 
4U 0614+91 with M ~ 0.01M Edd (IKuulkers et al.ll2010l) . 
and 4U 1608-522 where the accretion rate at the time 
of the superburst was high, but where the average rate 
over the previous years was M ~ O.OlMEdd (jKeek et al.l 
2008). The a-parameter, the ratio of the accretion flu- 
ence between normal burst s and the fluence of a burst, 
is typically high: a ~ 1000 (jin 't Zand et al.ll2003h . This 
indicates a relatively large part of the accreted material is 
burned in a stable manner instead of in bursts, and this 
may be necessary to achieve high enough carbon frac- 
tions. No superbursts have been observed from sources 
that only have stable burning and no bursts, though 
lower limits on the possible r ecurrence time have been 
determined (|Keek et alJ 120061 ). Although bursts reduce 
the carbon content of the envelope in the production of 
heavy elements, it has been suggested that the heavy ele- 
ments are necessary for reducing the thermal conductiv- 
ity, insuring that the superburst ig nition is reached at the 
obser ved depth in the envelope ([Cumming fc Bildstenl 
120011) . 



2 



From fits of superburst-decay mode ls to observed light 
curves ([Cumming fc Macbeth! I2004T ) . iCumming et al.l 
( 2006) deduce that superbursts ignite at a column depth 
of y ~ 10 11 — 10 12 gcm -2 in a layer with a carbon 
mass fraction of X\2 — 20%. It is a challenge for mod- 
els to explain these ignition column depths. The car- 
bon mass fractions are higher than what one-dimensional 
models that include large nuclear netw orks predict to 
be present in the ashes o f normal bursts (Woosle v et al.1 
[200l Ifisker et alj l200l . iCooper et all (|2006f ) suggest 
that the companion stars of superbursters donate mate- 
rial with a CNO content that is four times higher than 
solar. 

Superbursts ignite close to the outer crust, and as such 
are sensitive to the thermal p roperties of th e crust, which 
are not yet well understood (jBrownl [20041) . In turn, the 
temperature of the crust depends on neutrino-cooling 
in the neutron star core, which is also ill-constrained. 
Therefore, superbursts provide an observational measure 
of the thermal properties of the outer cru st, and constrain 
the physics in the crust and the core (jCumming et all 
2006; IPaee & Cum"mTnl[2005l) . 

The start of the superburst was observed only in eight 
cases. In six of these, a short precursor burst is de- 
tected. For the other two superbursts, the data was 
not of sufficient quality to exclude the presence of a 
precursor, with the possible exception of 4U 1608-522, 
although the detection of the supe rburst onset mus t 
be regarded tentative for th is source (IKeek et al.l [20081. 
iWeinberg fc Bildstenl (|2007l) explain the precursor as the 
result of a shock generated by the superburst ignition. 
This shock travels outwards through the envelope and 
triggers the ignition of either a helium-rich layer or an- 
other carbon-rich layer. The resulting flash is observable 
as the precursor burst. 

In this paper we create a series of one-dimensional 
models of the neutron star envelope, where for the first 
time we self-consistently build up a carbon-rich layer at 
rates similar to the observed accretion rates. We fol- 
low the carbon burning during several consecutive su- 
perbursts. The dependence of observable properties of 
the bursts on crustal heating is investigated. A possible 
hydrogen or helium-rich atmosphere is not modeled in 
this paper. 

2. NEUTRON STAR ENVELOPE MODEL 
2.1. Stellar evolution code 

We employ the one-dimensi onal hydrodynamic s stel- 
lar evolution code KEPLER (|Weaver et all 11978! ) . We 
use a version of KEPLER t hat differs from the version 
used in recent st udies (e.g., iWooslev et al.l 12003 . 120041 : 
iHeger et al.ll2007h in the accretion scheme and the opac- 
ities that are used. We model the neutron star enve- 
lope on a one-dimensional Lagrangian grid in the radial 
direction, under the assumption of spherical symmetry. 
The grid points represent the boundaries between con- 
centric shells, that each have a certain mass, chemical 
composition, temperature, density, luminosity and ra- 
dial velocity. Alternatively, a model could be considered 
a local 'wedge' of the neutron star, that would be well- 
approximated in a plane-parallel geometry. 

Zones are added and removed in order to maintain an 
optimal grid for resolving gradients in all quantities, such 



that temperature, density, and radius vary from one zone 
to the next by at least 10% and at most 25%. Further- 
more, zones are not removed if they extend over 0.02 in 
logy, where y is the column depth. The effects of differ- 
ent rezoning criteria was tested in a limited number of 
models; the most important calculated properties such as 
burst recurrence times and energetics varied by at most 
a few percent. The mass of each zone as well as the size 
of each time step are recorded, such that small values are 
not lost due to numerical precision. 

We implicitly solve the eq uations of mass, e nergy, 
and momentum conservation ([Weaver et al.l lT978). The 
equation of state allows for (non-) degenerate and (non- 
)relativistic electrons. 

To follow the chemical evolution we have the use of 
two networks of nuclear reactions. An adaptive network 
follows a large number of rea ctions among hundreds of 
isotopes (jWoosley et alJ [2004( 1. Because this network is 
computationally expensive, most of our calculations only 
use an approximatio n network consisting of 19 isotopes 
(| Weaver et al.l [l978h . It includes the carbon fusion re- 
actions as well as photodisintegration. Comparison of 
superburst models created using either network shows a 
3.3% shorter recurrence time for the model with the ap- 
proximation network, and a 3.1% lower burst fluence. 
This indicates that the large network generates 0.2% 
more energy per unit mass than the approximation net- 
work. There is no notable difference in the li ght curves. 

We take into account n eutrino ene rgy loss ([Itoh et "all 
I1996T ). ra diative opacity l|Ibenlll975[ ). and electron con- 
ductivity (jltoh et al.ll2008t) . 

We consider convection using the Ledoux criterion, as 
well as semiconv ection and thermohaline mixing (e.g., 
IHeger et all [2000). The induced mixing of the chemical 
composition is implemente d as a diffusiv e process using 
mixing-length theory (e.g.. lClavtonlll968j ). Rotation and 
magnetic fields are not considered in these models. 

2.2. Accretion and decretion 

Previous studies of X-ray bursts with the KEPLER 
code implemented accretion by increasing the pressure 
at the outer zone over time to simulate the build-up 
of a column of m aterial (e.g., IWooslev fc Weaver! 119841 ; 
iTaam et alj 119961 ) . When this pressure reached a cer- 
tain value, an extra zone containing the accreted mass 
was added on top of the model. This induced a momen- 
tary reduction of the time step as well as an artificial dip 
in the light curve. In the present study we employ an 
improved accretion scheme that solves these issues, al- 
lowing for larger time steps between subsequent bursts, 
and producing light curves without the aforementioned 
artifacts. 

Mass accretion is implemented by increasing the mass 
of one zone in the model at each time step at the mass ac- 
cretion rate. The zone is selected at a pre-defined column 
depth such that it lies above the region where thermonu- 
clear burning takes place, but far enough below the sur- 
face that the mass added to it is small compared to the 
layer above (this avoids constant rezoning of the small 
surface zones). Once the mass of the zone reaches a cer- 
tain limit, it is split together with one neighboring zone 
into three zones, conserving energy, momentum, compo- 
sition, and gradients. The chemical composition of the 
zone and all zones above, up to the surface, is advected 
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to account for the composition of the accreted material. 
Furthermore, the radial positions of the zones above the 
mass addition point are adjusted, and the energy gained 
from compressional heating of the accreted material is 
taken into account. 

Increasing the mass of a model leads to increased neu- 
trino emission near the bottom. To avoid this we main- 
tain a constant total mass for the model, by decreasing 
the mass of the inner zone at the same rate as at which 
mass is accreted. The radius of the inner boundary is 
kept fixed, and all other zones are moved downward, such 
that the density in the inner zone is conserved. Once the 
first zone's mass is reduced below a certain limit, the 
three inner zones are merged into two, again conserving 
energy, momentum, composition, and gradients. 

2.3. Substrate 

The inner part of the models is formed by an iron sub- 
strate, on top of which the carbon-rich superburst fuel 
is accreted. Heat generated in a burst can diffuse into 
the substrate, and be released toward the surface on a 
longer time scale. This ensures a correct long-term light 
curve. We performed tests that show that the substrate 
should contain at least an order of magnitude more mass 
than the burst ignition column. At low accretion rates 
the ignition column depth is relatively large, requiring 
the substrate to be located deeper. 

The substrate lies below the superburst ignition depth, 
and reaches into the outer crust, where neutrino emission 
becomes increasingly important. If we choose too thick a 
substrate, most of the luminosity at the inner boundary 
will be dissipated as neutrinos. This is especially a prob- 
lem for models with high accretion rates, which have a 
relatively high crustal heating and thence larger neutrino 
losses. 

The wide range of ignition column depths and amounts 
of neutrino losses pose constraints on the substrate mass 
that vary as a function of the mass accretion rate. At 
the lowest rates we choose the substrate to have a mass 
of 2 ■ 10 28 g, and 2 • 10 26 g at the highest rates. As a test, 
we create several models with the same accretion rate 
and varying substrate sizes. The changes in the burst 
parameters such as the recurrence times is at most a few 
percent for the selected substrate sizes. 

2.4. Crustal heating 

The amount of crustal heating of the envelope depends 
on the nuclear reactions in the neutron star crust, the 
crust's thermal conductivity, and on the neutrino cooling 
in that layer and in the core. The processes in the crust 
are not calculated explicitely, but the resulting heating 
of the envelope is emulated by a fixed luminosity at the 
inner boundary. For each model, we assume a certain 
heat flux per accreted nucleon Qt>- Combined with the 
accretion rate, it specifies this luminosity. 

The inner part of our models, the substrate, reaches 
into the crust, and the luminosity that reaches the super- 
bursting region is reduced by neutrino emission. Because 
we wish our results to be independent of our prescription 
of crustal neutrino cooling, we report an effective Qb, 
that is corrected for neutrino emission in the substrate. 

2.5. Relativistic corrections 



The code we employ uses Newtonian gravity (calcu- 
lated for each zone), whereas for neutron stars general 
relativistic (GR) corrections are significant. To take 
these corrections into account, we can state that our re- 
sults are applicable to any combination of neutron star 
mass and radius that give a GR gravitational acceler- 
ation equal to the Newtonian acceleration employed by 
the code. The full details of the GR corrections are avail- 
able in Appendix [B] Here we give one example, but note 
that the results of the models are valid for any combina- 
tion of mass and radius that yields the same value of the 
gravitational acceleration as used in this study. 

An input mass of 1.4 M@ and radius of 10 km yield 
a local Newtonian gravitational acceleration throughout 
the envelope of g ~ 1.87 • 10 14 cms~ 2 . Using the same 
mass and a larger radius of 11.2 km, one obtains the same 
value of the gravitational acceleration, but now including 
GR corrections. So the results of the Newtonian model 
are valid for a GR model with increased radius. Because 
of the larger radius, the luminosity from our model has 
to be increased as well, by a factor 1.12 2 ~ 1.26. This 
mass and radius imply for an observer at infinity a gravi- 
tational rcdshift of 1 + z ~ 1.26. The observed luminosity 
is reduced by a factor 1.26, and the observed ratio of the 
accretion luminosity and the Eddington limit is scaled by 
a factor 0.99. The GR corrected global mass accretion 
rate for an observer at inifinity is the same as the input 
(non-redshifted) accretion rate M. 

The results presented in this paper do not contain these 
corrections unless indicated otherwise (e.g., Sect. 13 .6j ). 

2.6. Light curves 

Light curves a re generated using the luminosity in the 
outer zone fe.g.. iTaam et al.1ll996| ). This zone typically 
extends orders of magnitude in column depth deeper be- 
low the neutron star surface than the photosphere. The 
surface zone, therefore, has a much longer thermal time 
scale than the photosphere. For our models this typ- 
ically means that thermal diffusion cannot change the 
light curve faster than on a thermal time scale of ~ 10 -4 s 
for an outer zone of 10 16 g. Dynamic processes such as 
shocks, however, can heat the outer zone much faster, 
producing variations in the light curve on shorter time 
scales. We do not correct for the time it would take to 
transport heat through the outer zone to the 'real' pho- 
tosphere, which is a reasonable approximation for the 
dynamic processes because of the short spatial distance 
to the photosphere. 

As explained in the previous subsection, no GR cor- 
rections are applied to the light curves. 

2.7. Initial model setup 

For the inner boundary we set the radius to 10 km 
and the enclosed mass to 1.4 M® (using Newtonian grav- 
ity). The outer boundary is initially set at a column 
depth of y — 10 9 gcm~ 2 . Once accretion is turned on, 
zones are quickly added such that the boundary is at 
y = 10 3 gcm~ 2 , which corresponds to the outer zone 
having a mass of - 10 16 g. Note that we refrain from 
resolving the photosphere at y ~ lgcm~ 2 , because this 
would require very light zones that display unphysical 
behavior in the presence of shocks. 

The model initially consists of the iron substrate (50 
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Figure 1. Mass accretion rate M in units of the Eddington- 
limited rate, Mg^di an d effective crustal heating parameter, Qbi 
for 78 models. For each model we indicate whether carbon burn- 
ing proceeds in a stable manner or as a superburst. The lines 
correspond to differen t models for crustal heating (Fig. 18 from 
Gumming ct al. 2006): models with higher neutrino cooling have 
lower values of Qb (solid lines); a model with a highly impure crust 
(dot-dashed line). 

zones) and five zones containing a mixture of 80% 56 Fe 
and 20% 12 C, which is later used as accretion com- 
position. This is the typical mass fraction of carbon 
ICumming et al.l ((2006) found from fits to observed light 
curves of hydrogen accreting superbursters, and 56 Fe is 
the most abundant is otope in the ashes of hydrogen- 
rich X-ray bursts (e.g., Wooslcv ct al. 2004). Because of 
the different heat sources (crustal, compressional heat- 
ing) and sinks (surface radiation, neutrino emission), the 
model must be brought into thermal equilibrium before 
the simulation is started. The model is evolved by the 
code over a period that is much longer than the thermal 
time scale of any zone. In this period we do not consider 
nuclear burning and mixing processes. With respect to 
accretion, we do not change the mass and composition 
of the model, but we do advect compressional heating 
throughout all zones. Crustal heating is applied as well. 
Once the model is in thermal equilibrium, we reset the 
simulation time to zero, and enable accretion fully, as 
well as nuclear burning and mixing processes. Because 
of the accretion of mass, zones are added: a typical model 
contains approximately 400 zones, about 50 of which are 
located in the substrate, and another 50 form the outer 
region where accretion is implemented. 

We create models with different values for Qb and 
for the mass accretion rate M, expressed as a frac- 
tion of the Eddington limited rate AfEdd- In this pa- 
per we use the Eddington limit for an atmosphere of 
solar composition on a neutron star of 1.4M with a 
11.2 km radius, which corresponds to an Eddington lu- 
minosity of iEdd = 2.5 • 10 38 ergs _1 and accretion rate 
of MEdd = 1-96 ■ 10~ 8 Mq year -1 . Taking into account 
the gravitational redshift, the observed values at in- 
finity are Z/ E ddoo = 1.6 • 10 38 ergs -1 and M E ddoo - 
1.56 • lO- 8 M year -1 . 

3. RESULTS 

3.1. Stable/unstable ignition 

We create a series of models for different values 
of the mass accretion rate and crustal heating. We 
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Figure 2. Superburst ignition column depth j/i gn as a function of 
the crustal heating parameter Qb for series of models with a fixed 
mass accretion rate indicated as a fraction of M^dd- 
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Figure 3. Ignition (burning) temperature Ti gn , density pi gn and 
column depth j/j gn for models of superbursts (stable carbon burn- 
ing), in a wide range of mass accretion rates and crustal heating. 
The relation between pj gn and y\ &n used here (a power-law fit to 
the numerical results) is accurate up to 1.2%. 

vary the mass accretion rate in the range where su- 
perbursts are observed 0.0 1 < M/M Edd < 1.00 
(e.g. JKeek fc in 't Zan d 2008), and we vary the amount 
of crustal heating between the minimum and maxi- 
mum val ues suggested in the literature: 0.1 < Qb < 
1.5 (e.g..lHaenseYfe Zdunlkll2003l ICumming et al.ll2006l : 
IGupta et alJ I2007IK In F"igT[Tl we indicate for which val- 
ues of these parameters the models exhibit superbursts 
or stable carbon burning. We include the different pre- 
dictions for Qb as a function of M from ICumming et alJ 
(2006), with the exception of the model that includes 
Cooper pairs, as its neutrino emissivit y was shown to be 
overestimated (|Leinson & Perez 2006). 

Ignition occurs in our models at a column depth of 
?/ig n ~ 10 10 gcm~ 2 for the models with the highest accre- 
tion rate and crustal heating, and at yi gn ~ 7T0 13 g cm -2 
for the coolest models with the lowest accretion rate 
(Fig. [21 O- We determine y; gn in our models from the 
location of the peak temperature just before the start of 
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Figure 4. Superburst recurrence time, t TCCUT , as a function of the 
crustal heating parameter, Qb, for series of models with a fixed 
mass accretion rate, indicated as a fraction of M E dd- The dotted 
extrapolations continue after the hottest bursting models to the 
Qb value where the first stable model was found in our grid. For 
three sources we indicate the shortest observed recurrence time (at 
arbitrary Qb)- The average lower limit from Be ppoSAX WFC da ta 
for nine sources is indicated by a vertical line (Keek et al. 2006). 
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Figure 5. Energy generation rate, E U uc, surface luminosity, L, 
and neutrino luminosity, L„, as a function of time from the start 
of one superburst up to the onset of the next, for a model with 
M = 0.30 M Edd and Q b = 0.13 MeV nucleon -1 . 

the runaway (in case of unstable burning), or the peak 
energy generation rate (in case of steady-state burning), 
which in both cases identifies the bottom of the carbon- 
rich layer. 

The bursting models exhibit recurrence times of several 
days up to thousands of years (Fig. |4])- A given recur- 
rence time can be reproduced by a relatively hot model 
with a certain accretion rate or a colder model with a 
somewhat higher accretion rate. 

3.2. Thermonuclear burning 
To illustrate the thermonuclear burning processes dur- 
ing a superburst we consider a model with M = 
0.30M E dd and Q b = 0.13 MeV nucleon" 1 . The energy 
generation rate is highest ~ 2 • lCP 7 s after the ther- 
monuclear runaway, and decreases roughly as i -1 over 
the course of ~ 10 6 s (Fig. [5]). Afterwards, the energy 
generation rate rises again due to increased carbon burn- 
ing in a newly accreted fuel layer, leading up to the next 
superburst. 
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Figure 6. Total mass for a selection of isotopes as a function of 
time from the start of one superburst up to the onset of the next, 
for a model with M = 0.30 M Edd and Q b = 0.13 Me V nucleon - 1 . 
For Fe the mass at t = is subtracted: only the produced mass 
is shown. The Fe mass decreases toward the end because of mass 
removal at the inner zone of the model (Sect. 12.21 . not because of 
nuclear reactions. 
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Figure 7. Mass fraction per isotope as a function of column depth, 
y, approximately 10 4 s after the start of a superburst for a model 
with M = 0.30 M Edd and Q b = 0.13 MeV nucleon -1 . The region 
below y = 10 7 gcm -2 is omitted, as it has constant carbon and 
iron fractions. The column depth at which the superburst ignited, 
j/ig n , is indicated by the vertical dotted line. 

At the superburst onset a large fraction of the avail- 
able carbon burns through the 1 C( 12 C,a) 20 Ne reaction 
(Fig. [5]). Subsequent a-capture reactions produce heavier 
isotopes, such as 24 Mg, 28 Si, and 32 S. Photodisintegra- 
tion causes the release of more a-particles, whose cap- 
tures create iron-group elements (e.g., iron and nickel). 
Electron-captures onto nickel produce iron, which is the 
most abundant element in the superburst ashes after sev- 
eral seconds. Note that in the employed approximation 
network the only iron isotope is Fe, whereas calcula- 
tions with a large network confirm that 56 Fe is the most 
abundant isotope. 

Approximately 10 4 s after the burst start, the total 
amount of carbon increases again when accretion adds 
carbon faster than residual burning can take it away 
(Fig.[6|). At this time there is a layer of pure iron directly 
above the ignition depth, that accounts for over 90% of 
the mass of the superbursting layer (Fig. [?])■ This is the 
composition that forms the outer crust. In the outer part 
of the envelope, photodisintegration was less efficient due 
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Figure 8. Ignition column depth, yi gn , and the total burst energy, 
Eburat, emitted at the surface for all models. The linear relation 
below j/i gn < 10 12 gem -2 is indicated by the dotted line. 

to the lower temperature, resulting in ashes that are more 
rich in 20 Ne, 24 Mg, 28 Si, and 32 S. The new fuel piles on 
top of that layer. During the next superburst these iso- 
topes burn to iron group elements. Note that models 
with stable carbon burning do not undergo photodisin- 
tegration, and there the outer crust will be enriched with 
20 Ne, 24 Mg, 28 Si, and 32 S. 

Part of the generated energy leaves the envelope as 
photons from the surface, and part is lost in neutrinos 
(Fig. [5]). Neutrino emission is strongest at larger col- 
umn depths in the substrate. The total burst energy 
emitted in photons at the surface, -Ehurst, for all burst- 
ing models follows a linear relation for ignition below 
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10 12 gcm~ 2 (Fig. [SJ. At larger depths E huTSt 



drops below this relation, because an increasing part of 
the energy is emitted as neutrinos from the substrate 
(the crust). The maximum -Eburst in these models is 
1.1 • 10 43 erg. 

3.3. Shock and mixing 

To study the hydrodynamic processes during a super- 
burst we again consider the model with M = 0.30 M^dd 
and Qb = 0.13 MeVnucleon -1 . When the thermonuclear 
runaway occurs at the bottom of the carbon-rich layer, 
the burning initially proceeds as a detonation. A com- 
bustion wave moves outward, creating a shock. After 
several microseconds the combustion wave slows down, 
and burning spreads as a deflagration. The shock con- 
tinues to travel toward the surface on a microsecond time 
scale (Fig. [5] top). The top layers are pushed outward, 
and subsequently fall back on a dynamical time scale of 
approximately 10 -5 s. Afterwards the surface undergoes 
a damped oscillation. In Fig. [5] we only show the enve- 
lope down to y ~ 10 s gem -2 , which provides the most 
insight into the dynamic processes (see also Sect. 13. 51) . 

As the outer layers fall back, most of the kinetic en- 
ergy is dissipated into heat at a depth of y ~ 10 8 to 
10 10 gcm -2 . Heating by the shock and the fall-back in- 
duces some carbon burning in this region, leading to two 
regions of thermonuclear burning (Fig. [5]top). 

Several 100 seconds before and after the superburst 
onset, convection mixes the composition in the envelope 
(Fig. [H] middle). Briefly, at the thermonuclear runaway, 
the convective region reaches close to the surface. Af- 
ter the burst, burning continues at a very low rate at 
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Figure 9. On three time scales: energy generation/loss (color 
scale) in the neutron star envelope as a function of time in a 
short interval around the superburst onset for a model with M = 
0.30 AfEdd and = 0.13 MeV nucleon - 1 . Green hatching indi- 
cates convection, red cross hatching semiconvection. 

the bottom of the freshly accreted layer (Fig. |9] bottom). 
During this time no convective mixing takes place. Once 
the ignition column depth is reached, the next superburst 
occurs. 

The compositional gradient created by the superburst 
induces thermohaline mixing in the envelope. This mix- 
ing is, however, many orders of magnitude slower than 
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that due to convection at the time of burst onset. 

3.4. Light curve 

We ge nerate light curves for all bursting models 
(Fig. I A 1 p . We compare a selection of light curves by 
taking the coldest model (Qb — 0.15 MeV nucleon~ 1 ) in 
a wide range mass accretion rates (Fig. [TO)) . The curves 
consist of several components: a shock breakout peak, a 
precursor, a transition to the superburst peak, followed 
by a two-part power law decay. Not every model exhibits 
each component. After the superburst peak, the decay 
proceeds as t -021 . Following the break, the decay steep- 
ens to i~ 136 . The models with the lowest mass accretion 
rates, which are the coldest models with the largest igni- 
tion column depth, have the longest decays. The t~ 21 
power law forms an upper bound to the light curve. To- 
ward higher accretion rates, the time spent in the £~ ' 21 
part is smaller, until it becomes absent at the highest 
rates, and there is a direct transition from the peak to 
the t~ 136 decay. 

Whereas the models with lower accretion rate exhibit 
a precursor, at the highest accretion rates — the hottest 
models — it is absent. Colder models have longer precur- 
sor bursts of up to ~ 10 2 s. We find precursors as short 
as ~ 10 _1 s. All precursors reach the Eddington limit 
and cause radius expansion. Depending on the duration 
of the precursor, the transition to the superburst 'peak' 
around ~ 10 3 s can exhibit a drop in luminosity below 
the peak value. The models with longer precursors lack 
this dip. 

To study the effect of crustal heating on the light curve, 
we compare a series of light curves of simulations with 
the same mass accretion rate, M = 0.30MEdd> but with 
increasing Qb (Fig. fTTj) . For increasing Q_, the duration 
of the superbursts increases from 2.0 • 10 5 s to 6.7 • 10 5 s, 
and the time spent in the i~ 21 phase reduces. The 
hottest models again lack the precursor. Comparing the 
models with Qb = 0.13 and Qh = 0.17, the latter has 
a shorter duration precursor, as well as a deeper dip in 
the light curve between the precursor and the main peak 
after 1000 s. 

Some light curves show a shorter precursor phase 
than expected, exhibiting instead a small bump im- 
mediately after the precursor (e.g., some models with 

M = 0.10 M E dd in Fig. IXljl . This is due to the relax- 
ation of the outer atmosphere following the end of the 
radius expansion phase, and may be attributed to poor 
resolution at the surface of the models. It is not caused 
by burning or mixing processes. 

We define the exponential decay time, t cxp , of our 
light curves as the time it takes from the superburst 
peak around t = 10 3 s to drop one e-fold in luminosity 
(Fig. [T2")) . It ranges from 18 minutes to 5.2 hours. 

3.5. Precursor burst 

As noted bv lWeinberg fc Bildstenl (|2007r i , the overpres- 
sure of the shock is larger at lower column depths. Hence, 
at lower column depths the shock induces more heating 
and a larger radius expansion. To illustrate this, we com- 
pare the light curves of two models that have both M = 
0.30 AfEdd and Qb = 0.13MeVnucleon _1 , but one model 
extends to a column depth of y ~ 10 s gem -2 , and the 
other to y ~ 10 3 gcm~ 2 (Fig. [l~3|) . The model with the 



more extended envelope has a shock breakout peak that 
has an approximately 50 times faster rise and a super- 
Eddington luminosity of 0.89 • 10 41 ergs _1 , whereas the 
other model's shock peak reaches only 1.4 • 10 37 ergs _1 . 
After the shock breakout, the model with the more ex- 
tended envelope has stronger radius expansion and larger 
variations in luminosity. The latter reaches a 1.7 times 
higher value than for the other model. The stronger 
shock heating leads to a lower opacity, which increases 
the Eddington limit, allowing for higher surface lumi- 
nosities. 4 seconds after the start of the superburst, the 
radius expansion phase ends. The more extended atmo- 
sphere model displays a sharper drop in luminosity. The 
smoother luminosity decrease of the other model may 
lead to the interpretation that the precursor in this case 
has a duration that is several seconds shorter. After the 
precursor, the light curves of the two models are virtually 
identical. Therefore, both the duration and the peak lu- 
minosity of the shock breakout and the precursor depend 
greatly on the extend of the atmosphere. Note that the 
neutron star photosphere is expected at a column depth 
of approximately y ~ lgcm -2 , and it is likely that at 
that column depth the shock breakout peak as well as 
the precursor properties are different from our results. 

The shock heating of the outer layers induces some car- 
bon burning, but the amount of heat generated by the 
nuclear reactions is too small to substantially alter the 
light curve, as we checked by comparing to a model where 
burning was disabled in that region. Therefore, the pre- 
cursor is in these models virtually completely powered 
by shock heating. 

The hottest models at the highest accretion rates do 
not show precursors at all (Sect. I3.4[) . We compare 
the model with a precursor from Fig. i to a hotter 
one (Qb = 0.22 MeVnucleon -1 ) that lacks a precursor 
(Fig. [T4|) . The hotter model's burst has a shallower ig- 
nition depth, and is, consequently, less powerful. The 
shock causes only minimal radius expansion, and does 
not provide enough heating to produce a precursor burst 
or to ignite carbon close to the surface. In this case there 
is only one region of carbon burning (Fig. I14[) . 

3.6. Comparison to 4U 1636-53 

We compare several of our simulated light curves to the 
observed light curve from the superburst of 4U 1636-53, 
that was observed with the Proportional Counter Ar- 
ray (PCA) on-boa rd the Rossi X-ray Timing Explorer 
(R XTE1 (Fig [TBI i Strohmaver fc Markwardi) 1200a see 
also|Kuulkers 2004; IKuulkers et al.ll2004f) . The PCA con- 
sists of five Proportional Counter Units (PCUs) . We use 
standard 1 data from PCU 2 in the full instrument band- 
pass, and correct it for dead time following the prescrip- 
tion from the RXTE Cookboo __ We subtract the persis- 
tent emission as measured from the end of the last orbit 
that we consider. In this we assume the persistent flux to 
remain constant during the superburst, which is proba- 
bly not the case: in the day preceding the superburst, the 
persistent flux varied by around 10 2 counts s^ 1 PCU -1 . 
Because the superburst lasts longer than the orbit of 
RXTE, the observation was interrupted three times by 
Earth occultations, resulting in gaps in the light curve. 

1 http://heasarc.gsfc.nasa.gov/docs/xte/recipes/pca_deadtime.html 
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Figure 10. Light curves of superburst models with Qj, ~ 0.15 MeV nucleon for a wide range of mass accretion rates. 
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Figure 13. Light curves of superburst models with M = 
0.30AfEdd and Qj, = 0.13 MeV nucleon" 1 with different column 
depths for the outer zone. 
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Figure 12. Exponential decay time, t eX p, vs. crustal heating, Qbi 
for series of models with a certain mass accretion rate, indicated 
as a fraction of MEdd- 

The simulated light curves are constructed taking into 
account the PCA's instrument response and astrophys- 
ical effects. We use the surface radius and tempera- 
ture from our models to calculate the blackbody emis- 
sion from the superbursts. The temperature is increased 
by a typical color correction factor of 1.5 to account for 
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Figure 14. For a model without a precursor (M = 0.30AfEdd 
and = 0.22 MeV nucleon -1 ) the energy generation/loss (color 
scale) in the neutron star envelope as a function of time in a short 
interval around the superburst onset. Green hatching indicates 
convection. 
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53 to light curves of superburst models with M = 0.10 and 
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deviation from a pure blackbody spectrum due to Comp- 
ton scattering close to the neutron star surface (e.g., 
iSuleimanov et all 120111 ). We apply a gravitational red- 
shift of 1+z = 1.26 for a neutron star with a gravitational 
mass of 1.4 M Q and a radius of 11.2 km in the local rest 
frame (see Sect. 12.51 Appendix [B|). The effect of interstel- 
lar absorptio n by hydrogen is taken into acc ount using 
the model bv lMorrison fe McCammon (1983), that uses 
solar abundances from I Anders fe Ebiharal (jl98"2l) , using 
a hydrogen column of 2.5 • 10 21 cm -2 (jAsai et al.l 120001) . 
We take into account the effective area of the PCUs at 
different energies using the table provided with the soft- 



ware package PIMMS version 4.2. The curves are scaled 
such that the superburst peak fluxes match at t ~ 800 s. 

We obtain a measure of the persistent luminosity dur- 
ing the month preceding and following the superburst 
from flux measurements obtained with the PCA on 
RXTE and the Wide-Field Cameras (WFCs) on-board 
BeppoSAX, that are collected i n the Multi-Instru ment 
Burst Archive (MINBAR; e.g., iKeek et alj 120 lOl) . A 
bolometric correction is available for several orbits. We 
use the mean value: 1.4. By comparing to the Eddington 
luminosity for a 1.4 M@ neutron star with a 11.2 km ra- 
dius and an atmosphere of solar composition (Sect. l2~7|) . 
we find that the accretion rate was 0.12 MEdd with a root 
mean-squared deviation of 10%. 

For models with an accretion rate of 0.10 MEdd, the 
scaling factor for our simulated curves is approximately 
0.42. This factor can be explained by a 54% larger dis- 
tance to the source than the 6 kpc that we assumed 
([Galloway et al.l 120081 ). These models predict a much 
longer decay time than observed (Fig. ITS"]) . The best 
fit is provided by a model at higher accretion rate: 
Q b = O^lMeVnucleon- 1 and M = 0.40 M Ed d- For 
those models the scaling factor is approximately 0.63, 
which can be explained by a 26% larger distance to the 
source. The best fit model lacks a precursor burst. 

4. DISCUSSION 
4.1. Superburst models 

We create 86 models of a neutron star envelope with 
accretion of carbon-rich material in the range of observed 
accretion rates, assuming a range of values for the crustal 
heating parameter Qb (Fig.Q]). We follow the thermonu- 
clear burning of the accreted carbon, which proceeds as 
flashes (superbursts) in some cases, whereas at higher 
values of Qb burning becomes stable. We compare the 
amount of crustal heating of our bursting models to mod- 
els of Qb as a function of M jC ummmg eTaT1[20(il . Only 
the lowest curve that spans the entire range of mass ac- 
cretion rates lies within the range of Qb where we find 
unstable carbon burning. This implies that a high neu- 
trino emissivity of the neutron star core is favored. 

Note that superbursting sources mostly accrete hydro- 
gen or helium-rich material, which create carbon-rich 
ashes from thermonuclear burning. By directly accret- 
ing the latter composition, we skip the computationally 
expensive hydrogen/helium burning, making it possible 
to simulate the long superburst recurrence times. Hy- 
drogen/helium burning may increase the temperature of 
the envelope, which can be modelled by an extra contri- 
bution to Qb- 

Carbon ignition occurs in our models at a col- 
umn depth yi gn between approximately 10 10 gcm~ 2 and 
10 14 gcm -2 (Fig. [21 Fig. [3]). The stable burning models 
extend to lower z/i gn than the bursting models, because 
at a given accretion rate stable burning requires a higher 
crustal heat flux than unstable burning, which leads to 
shallower ignition. 

In the relation between the ignition depth (or den- 
sity) and temperature (Fig. [3]), there is a down turn in 
the trend for yi gn > 10 13 gcm~ 2 , which is due to in- 
creased screening of the Coulomb barrier of the carbon 
ions at temperatures below T < 3 • 10 8 K and densities 
in excess of p > 10 9 gcm~ 3 (Salpete r fc van Horn! 11969} . 
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Figure 16. Light curve of the model with M = 0.30 M^dd antl 
Qb = 0.13 MeV nucleon - 1 , where the outer zone is at a column 
depth of ~ 10 s g cm" 2 . We indicate the different phases of the 
superburst. Several followup superbursts are shown. 

Yakovlev etlll l200rl see also, e.g., iBrown fc Bildster] 
1998). This is the transition from thermonuclear burning 
to pycnonuclear burning, which sets in at temperatures 
T < 10 s K. The ignition conditions are more uncertain 
in this regime, because the possible formation of a crys- 
tal lattice may requir e higher densities for carbon fusion 
(|Yakovlev et al.ll2006l ). 

When the ignition column depth exceeds y- lgn > 4 • 
10 11 g cm -2 we find superbursts to be powerful enough 
to dr i ve a shock to the surface (s ee also [Weinber g et al.l 
12001 IWemberg fc Bildstenl [2007h . At the start of the 
thermonuclear runaway, the burning time scale at the ig- 
nition depth becomes shorter than the dynamical time 
scale. A combustion wave moves outward, depleting the 
inner zones of carbon, while initiating a shock. This 
detonation phase lasts only a few microseconds, until 
the velocity of the combustion wave is sufficiently re- 
duced, and burning continues to spread as a deflagration. 
The shock then no longer follows the burning front, but 
speeds ahead toward the surface. 

4.2. Decay profile 

The superburst light curves consist of a shock break- 
out, a precursor, a rise to the 'superburst peak', and 
two power-law decay phases with L oc i~ - 21 and L oc 
£-i-36 (Fi g- [X^|) , These pow er law indices are consi stent 
wit h the results obtained bvlCumming et all (120061) and 
bv IWeinberg fc Bildstenl (|2007l ). ICumming fe Macbeth! 
(|2004j ) explain that the first power law component is 
due to radiative cooling when a c ooling wave propagate s 
from the surface inwards (see also ICumming et al.| [2006). 
Once this wave reaches the bottom of the carbon burning 
layer, the cooling transitions to the steeper power law. 
Hotter models with high accretion rates burn shallower 
columns, such that the L oc t~°- 21 phase lasts shorter in 
these models. For models with the highest values for M 
and Qh, this phase is absent in the light curve. After 
the superburst peak an immediate transition is made to 
L oc i- 136 . 

The L oc i -1 - 36 ends when burning of newly accreted 
carbon starts dominating the light curve, up to the ther- 
monuclear runaway of the next superburst (Fig. 1161) . 

Observationally, the duration of the decay is often 
measured by fitting an exponential to the light curve. 



Even though theoretically we expect the decay to fol- 
low a power law, the quality of the data is often such 
that an exponential still provides a good fit. Decay 
times have been observed between 0.7 and 6.0 hours (e.g., 
IKeek fc in 't Zandll2008| ). We determine the exponential 
decay time t exp for our light curves as the time it takes 
for the luminosity to drop from the peak by one e-fold 
(Fig. [T2). If this time falls within the L oc i~°- 21 phase, 
where all models share the same curve, i exp is mainly 
dependent on the peak luminosity: a higher peak means 
e^Lpeak is reached earlier on the L oc t~ Q 21 curve. This 
yields smaller values of i e xp for colder models, because 
they have larger peak luminosities. If t exp reaches into 
the L oc t~ 136 phase, t exp depends less on the peak, 
and more on the actual decay of the light curve. Colder 
models with lower accretion rates have longer i oxp in this 
case. The difference in behavior during the two power 
law decay phases makes it difficult to use the published 
observed i G xp to constrain model parameters. It is, there- 
fore, more instructive to fit a two-comp onent power law 
decay, if the data qua l ity allows for it (ICumming et al.l 
[2006tlKeek etHimM iKuulkers et alJl201ft Fig. 115)). 

4.3. Precursor bursts 

At the start of the thermonuclear runaway a shock 
travels from the ignition depth to the surface. The shock 
breakout is visible in the light curve as a short bright 
peak in the light curve, when during a fraction of a 
microsecond super-Eddington luminosities are reached 
(e.g., Fig. [131). Tne atmosphc re is pushed upwards, be- 
fore it falls back on a dynamical time scale of the or- 
der of 10 _5 s, and undergoes damped oscillations. The 
shock and the fall-back deposit heat for the different 
models at column depths between y ~ 10 8 gcm~ 2 and 
y ~ 10 10 gcm -2 . The cooling of these layers on ther- 
mal time scales of up to approximately 10 2 s is visible as 
a precursor burst. We find that the peak luminosity of 
the shock breakout and the amplitude of the oscillations 
depends strongly on the resolution of the atmosphere in 
the models. The shock over-pressure is greater closer to 
the surface ([Weinberg fc Bildstenl [2007h . where the den- 
sity is lower, causing a stronger shock breakout. Note 
that not all models produce a precursor: superbursts 
with yi gn < 4 • 10 n gcm~ 2 are not powerful enough to 
drive a shock. These models still reach a temperature 
of T ~ 5 • 10 s K at a column depth of y = 10 8 gcm~ 2 . 
This is too cold for carbon burning, but it is hot enough 
to tri gger the therm onuclear runaway of helium burning 
(e.g., Bildsten 1998), leading to a short hydrogen/helium 
precursor burst. 

Our models confirm the scenario of precursor s 
due to shock heat i ng (|Weinberg fc Bildstenl 120071 ). 
IWeinberg fc Bildstenl (2007) created models of super- 
bursts in a neutron star with a helium-rich envelope. 
They found that without the helium, carbon burning at 
lower depths produces a weak precursor. In contrast, 
we find that the shock heating is much stronger than 
any carbon burning at this depth, resulting in a bright 
precursor from heating alone. Their light curves exhibit 
a precursor that starts several seconds after the shock 
breakout. The shock leaves a the outer envelope isother- 
mal, delaying the precursor emissi on by a thermal time 
scale. IWeinberg fc Bildstenl (|2007f ) describe that at this 
point in the simulation they assume a new hydrostatic 
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equilibrium in the outer layers. Instead, we perform a 
fully self-consistent calculation that includes the fall-back 
of the outer layers, which disrupts the fiat temperature 
profile, and causes the light curve to change much faster 
on a dynamical time scale. Furthermore, most of the en- 
ergy of the shock went into the expansion of the outer 
layers, and is only converted into heat during the fall- 
back. This generates, therefore, not only a precursor 
more quickly on a dynamical time scale, but also a much 
more powerful precursor burst that can last up to 10 2 s. 

The light travel time around a neutron star with a 
10km radius is 1.1 • 10~ 4 s. Therefore, assuming instan- 
taneous ignition throughout the entire envelope, any de- 
tail in the light curve at shorter times, such as the shock 
breakout and the subsequent oscillations, will be smeared 
out. 

Only for seven superbursts the onset has been ob- 
served. Three times a precursor b urst was seen directly 
at the start of the superburst (IStrohmaver fc Brown 
2002t IStrohmaver fc Markwardtl l2002t lin 't Zand et al. 

2003) . For two superburst candidates from GX 17+2 
the onset was seen, but because of the q uality of the data 
only hints of precursors were observed (|in 't Zand et all 

2004) . The tentative observation of the superburst rise 
from 4U 1608-522 with HETE-2 did not exhibit a precur- 
sor. The data quality allows for only a weak precursor. 

The precursors of 4U 1636-53 and 4U 1254- 69 have a 
higher peak flux than the superburst. For 4U 1820-30 
the superburst itself was exceptionally bright, and the 
precursor's peak flux was approximately 15% lower. The 
light curves of the precursor bursts exhibit a double peak, 
similar to photospheric radius expansion (PRE) bursts. 
For none of these bursts, however, is spectroscopic infor- 
mation available at sufficient time resolution or of suffi- 
cient quality to confirm the PRE nature of these bursts. 
Nonetheless, the precursor bursts in our models all reach 
the Eddington limit, resulting in PRE. 

There is a tentative observation of the onset of 
the superburst from 4U 1608-522 with the WXM and 
FREGATE instru ments on-board the HETE-2 satellite 
(jKeek et a-l.l l2008h . At the end of an orbit there is an 
increase in the count rate visible for 50 s. No short pre- 
cursor can be discerned. The superburst occurred during 
a transient outburst, when the accretion rate exceeded 
0.10 MEdd- The time-averaged rate in the years before 
the superburst, however, was only 0.01 MEdd- Our mod- 
els at low accretion rates include precursors with dura- 
tions in excess of 50 s. This suggests the possibility that 
the entire flare observed with HETE-2 was part of the 
p recursor. 

iKuulkers etaU (|2002fl refer to a burst from KS 1731- 
260, that occurred 200 s before the superburst rise, as 
precursor. The 200 s is much shorter than the typical 
burst recurrence time for this source of several hours, as 
well as shorter than the superburst duration. The three 
other precursors, however, occurred immediately prior to 
the rise of the superburst. Also, this burst was relatively 
weak compared to most other bursts from this source, 
but several bursts with a similar peak flux have been 
observed. Therefore, we suggest that the burst preceding 
the superburst of KS 1731-260 is an ordinary burst that 
by chance was close to the superburst. Perhaps heating 
from the stable carbon burning before the superburst 



thermonuclear runaway caused the normal burst to ignite 
earlier, resulting in a relatively weak burst. 

4.4. Burning ashes and crustal composition 

Carbon burning and subsequent a-capt ure reactions 
create 20 Ne, 24 Mg, 28 Si, and 32 S (e.g., ISchatz et alJ 
I2003a|) . For stable burning models these isotopes com- 
bined with 56 Fe from the ashes of hy drogen /helium 
bursts (our accretion composition; e.g., IWooslev et al.l 
I2004D form the composition of the outer crust. In su- 
perbursting models, photodisintegration provides more 
a-particles for captures, creating iron group elements 
(Fig. [5]). Therefore, superbursting neutron stars have 
an outer crust composed of mostly iron. 

In the bursting models there is a region at a depth 
lower than the ignition column depth, where the temper- 
ature is insufficient for photodisintegration. Therefore, 
carbon burning in this layer creates 20 Ne, 24 Mg, 28 Si, 
and 32 S (Fig. [7|). Note that the most abundant element 
is iron from the accretion composition. The next super- 
burst ignites on top of this layer, and these isotopes burn 
to iron. A layer enriched in 20 Ne, 24 Mg, 28 Si, and 32 S in 
principle has a somewhat lower opacity than a pure iron 
layer, resulting in a larger ignition column depth yi gn . 
The layer composition is, however, dominated by iron, 
which limits the changes in opacity. We do not find a 
difference between y- lgn of the first superburst in a series 
and subsequent bursts that can be attributed to su ch a 
compositional inertia effect (cf. IWooslev et al.l [20041 ). 

4.5. Recurrence times 

The bursting models exhibit recurrence times of several 
days up to thousands of years (Fig. [4j) . We compare these 
to the three cases where more than one superburst has 
been observed from the same source. For these sources 
we derive the time averaged mass accretion rate from 
the persistent X-ra y luminosity as re ported in the MIN- 
BAR catalog fe.g lKeek et al If2010f ). expressed in units 
of the Eddington limited rate for our choice of neutron 
star parameters (Sect. 12. 7p . The typical uncertainty in 
measurements of the accretion rate are several tens of 
percents. Due to the presence of frequent data gaps, the 
observed re currence times ha ve to be regarded as upper 
limits (e.g., iKeek et"a"ni2010l for the case of bursts with 
short recurrence times). 

Ser X-l has a mean mass accretion rate of 0.31 MEdd, 
and superbursts obse r ved at least 2 .4 years apart 
(jCornelisse et al.1 120001 ; IKuulkers! l2009f ). The models 
with 0.30 MEdd have a maximum recurrence time of 
^rccur = 1.7 years (2.1 years for a gravitational redshift 
of z + 1 = 1-26), whereas models with 0.20 MEdd repro- 
duce the observed t rccur . Therefore, the accretion rate 
is somewhat lower than observationally inferred, or the 
actual i rccur of this source is shorter than the upper limit. 

Four superbursts ha ve been observed from GX 17+2; 
two only 15 days apart (|in 't Zand et a l. 2004). The mass 
accretion rate is unusually high for a bursting source: on 
average 1.21 MEdd- The models with a mass accretion 
rate of 1.00 MEdd exhibit recurrence times up to 11 days 
(13 days for z + 1 = 1.26). It is likely that models with 
a smaller accretion rate of ~ 0.9 MEdd will produce the 
observed recurrence time. This rate is within the uncer- 
tainty of the observed accretion rate. 
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The shortest time interval between t wo observed su- 
perbursts of 4U 1636-53 is 1.1 year (iWiinandsl [ 2001 1 : 
Strohmayer fc Markwardl l2002t iKuulkers et all 12004 : 
Kuulkers 2009). The time averaged mass accretion rate 

of 4U 1636-53 is 0.12 MEdd, but the models close to this 
rate produce a recurrence time of at least 7 years (9 years 
for z + 1 = 1.26). Part of this problem may be explained 
by the uncertainties in the measured mass accretion rate 
of several tens of percents. Other explanations for this 
discrepancy may be found in our assumptions for the 
carbon mass fraction and of the effective gravity in the 
neutron star envelope. 

For nine sources, most of which with accretion rates 
close to 0.10 MEdd, a lower limit to the superburst recur- 
rence time was derived from the BeppoS AX Wide-Field 
Camera (WFC) data (|Keek et al.ll2006[ ). The average 
lower limit of 60 days is indicated in Fig. 0] and is con- 
sistent with the model results for the sources with mass 
accretion rates up to 0.5A/Edd- 

The same recurrence time can be reproduced by models 
with a mass accretion rate that varies by several tens of 
percents, and different values of Qb- This spread in mass 
accretion rate is of the order of the uncertainty in the 
observed rate, which makes it difficult to constrain Qb 
from the observed recurrence times. 

4.6. Comparison to 4U 1636-53 

The most detailed superburst light curves have been 
observed with the PCA on RXT E: one from 4U 1636-53 
dStrohmaver fc Markwardtil2002fl and one from 4U 1820- 
30 (| Strohmayer fc Brown! 120021 ). The latter source is an 
ultra-compact binary system, which implies that the ac- 
creted material is hydrogen-deficient. The material that 
burns in the superburst is thought to be more carbon- 
rich than the composition that we assumed in our cal- 
culations. Its superburst is atypical, as the superburst 
peak reached the Eddington limit and exhibited radius 
expansion. For these reasons we do not compare to this 
source, and focus our attention on the superburst from 
4U 1636-53. 

The superburst from this source started with a short 
precursor that showed behavior consistent with photo- 
spheric radius expansion, and that decayed on a 10 s time 
scale. The superburst decay was observable for 5.5 hours. 
The average mass accretion rate of this source over the 
past decade is 0.12 MEdd fSect. l4~5l) . The models close to 
this rate exhibit a much longer decay time scale, even for 
high crustal heating (Fig. [TS] top). Also, the precursors 
of these models last longer than observed. The decay is 
best fit by a model with a four times higher accretion rate 
of M = 0.40 M Edd and with Q b = 0.21 MeVnucleon" 1 
(Fig. [T51 middle). This models, however, does not have 
a precursor, and the part leading up to the superburst 
peak is not well reproduced. This may be due to the 
fact that the atmosphere of 4U 1636-53's neutron star is 
probably hydrogen-rich, instead of carbon-rich as we as- 
sumed. While our modeled burst is not powerful enough 
to heat the atmosphere by a shock, the temperature in 
the atmosphere is high enough, T ~ 6 • 10 8 K, to ignite a 
hydrogen/helium precursor burst. 

If we were to assume that M = 0.40 M^dd is a good 
approximation, Fig. [T5] (bottom) illustrates how sensitive 
the decay depends on the crustal heating parameter Qb- 



The discrepancy between the observed light curve and 
the models at O.lOAfEdd adds to the problem we noted 
earlier that the recurrence times predicted by the models 
with 0.10 MEdd are substantially longer than the obser- 
vational upper limit (Fig. 3]). The answer to this problem 
may lie in the fact that we only considered one carbon 
fraction for the envelope, and one value for the effective 
gravity. Variations of these parameters could yield the 
shallower ignition implied by our models. Another possi- 
bility are multi-dimensional effects that we cannot model 
in our one-dimensional code. Pulsations at the spin fre- 
quency of the neutron star were observed during 800 s 
close to the superburst peak (|Strohmaver fc Markwardtl 
2002). This indicates that the emission was anisotropic 
during a quite long period after superburst ignition, 
which hints at the presence of multi-dimensional effects. 

5. CONCLUSIONS 

To study carbon flashes (superbursts) , we constructed 
86 one-dimensional multi-zone models of the envelope 
of a neutron star that accretes carbon-rich material, for 
different mass accretion rates and amounts of crustal 
heating. These are the first such models that were con- 
structed by following the accumulation of the fuel layer 
and the thermonuclear burning of carbon during a se- 
ries of superbursts. The stability of carbon burning is 
investigated as a function of the amount of crustal heat- 
ing. We reproduced the two- component power-law decay 
(jCumming fc Macbeth! [2004h . Not all models, however, 
exhibit the first component: the hotter models at higher 
mass accretion rates show a direct transition from the 
superburst peak to the second (steeper) power-law com- 
ponent. 

The superburst ashes that form the outer crust are 
primarily composed of iron. Carbon burning higher up 
in the envelope produces isotopes with mass numbers 
around 30. The next superburst ignites in this layer. We 
do not find a composit ional inertia effect, a s seen for hy- 
drogen/helium bursts (jWooslev et al.ll2004l) . because the 
layer composition is dominated by the accreted fraction 
of iron. In case of a larger carbon fraction of the fuel 
layer, however, such an effect may become important. 

We obtain a precursor burst due to shock heating, sim- 
ilar to lWeinberg fc Bildstenl (|2007l ). We find that heating 
by the shock and the fall-back of expanded layers is suf- 
ficient for a strong precursor, that starts approximately 
10~ 5 s after the shock breakout, instead of seconds. For 
hot models, at large accretion rates, the superburst is 
not powerfull enough to generate a shock and thus a pre- 
cursor from shock heating. At low accretion rates the 
precursors have durations as long as 10 2 s. This may ex- 
plain the lack of a short precursor in the observation of 
the onset of the superburst from 4U 1608-522. 

Comparing the model light curves to the superburst 
observations with the PCA on RXTE of 4U 1636-53, 
the models at the observationally inferred mass accretion 
rate overpredict the superburst duration and the recur- 
rence time. The best agreement is found with models 
at a three times higher accretion rate. The discrepancy 
may be caused by the values we assumed for the carbon 
fraction in the ocean and the effective gravity. This can 
be studied further with one-dimensional models. Alter- 
natively, it can be a sign of multi-dimensional effects, 
where a higher local accretion rate is responsible the ob- 
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served behavior. 

We studied the dependence of observables, such as the 
recurrence time and the shape of the light curve, on the 
amount of crustal heating Q^. While we show that these 
observables can depend strongly on Q^, the example of 
4U 1636-53 indicates that without a good agreement of 
the behavior as a function of mass accretion rate, it is 
difficult to constrain Qb- 
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APPENDIX 
A. MODEL LIGHT CURVES 

We present all light curves resulting from the models indicated in Fig. [1] (Fig. IA1[) . The curves do not contain 
corrections for the gravitational redshift near neutron stars. See Section [3.41 for further details. 

B. GENERAL RELATIVISTIC CORRECTIONS 
B.l. Introduction 

The thin surface layer of a neutron star in which the X-ray burst occurs can locally be treated in Newtonian physics 
for Ar <C r, where r is the radius of the star and Ar is the thickness of the surface layer. For neutron stars general 
relativistic (GR) effects are important, but the Kepler code, which is employed in this study, uses Newtonian gravity. 
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Fig. IA1I continued 

Therefore, we need to correct for GR effects to allow for proper interpretation and comparison with observations. 
This comprises two types of corrections: first, identifying the neutron star masses and radii that give same relativistic 
gravitational acceleration as the Newtonian acceleration that was used from the input mass and radius. Second, 
correct the results for time dilatation (of the light curve and decrease of the accretion rate; ~ (1 + z), where z is the 
gravitational redshift) and weakening of luminosity (~ (1 + z) 2 ). 

B.2. Translating Newtonian to GR 

When only Newtonian gravity is used in the calculation, neglecting the strengthening of gravity by a factor (1 + z), 
the result can still be interpreted as that of a star with a different mass, Ml, and radius, r 1; and correspondingly 
adjusted (smaller) redshift, zi, such that the GR acceleration equals the Newtion acceleration in the calculation. Below 
t he scaling laws for interpreting mass and radius are derived. 

iThornd (fl977T) gives the volume redshift factor 

V = l/y/l- 2GM/(c 2 r) , (Bl) 

which will be called 1 + z here. Relativistic gravitational acceleration is given by (e.g.. IWoodhousel 12007) 

g Ie l = -(! + z)GM /r 2 = GM / (r 2 ^/l - 2GM/(c 2 r)) . (B2) 

We now define a radius n and an actual gravitational mass, Mi, such that the GR gravity at this point equals the 
Newtonian gravity, g, at radius r, i.e., g r ci,i = 9'- 

GM jr 2 = GMi I (riVl - 2GMi/(c 2 n)) . (B3) 

This can be rewritten as 

ip 2 + 2ip(f - £ 4 = with ( = n/r, (f = Mi/M, ( = GM/(c 2 r) , (B4) 
where £ is the gravitational radius of the original problem. 

B.2.1. Given Mass 

Assuming a mass Mi = <pM, the physical solution of this 4 th order equation for £ is given by 



A= ^/2/9[B 2 /^ 2 -2^6) /(B( 2 ) , B= y9( 2 ^ + V3^ 3 ^16 + 27 C> 2 



(B5) 

The radius r\ at which the GR gravitational acceleration matches the Newtonian one for the assumed radius, r, is 
thus given by t\ = £r. The redshift factor z\ for radius r\ and mass Mi is given by 



1 



Using Equation (|B3|) one obtains 



+ zi = 1 / v /l-2GA/ 1 /( C 2n) = 1 jyf\ - 2^H . 



e/ v > = i+zi 



(B6) 
(B7) 
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Using these relations , the light curve for an observer at infinity has to be time dilated by 1 + z x . Due to the larger 
radius, the surface area is increased by £ 2 = ip(\-\-z\) and thus luminosity has to be scaled by £ 2 / (1+zi) 2 , i.e., decreased 
by a factor (1 + z x )/p. Similarly, the apparent accretion rate for an observer at infinity scales as £ 2 /(l + z x ) — <p, that 
is, does not need to be modified if M — M x . For a NS with 1.4 M© (gravitational) mass and 10 km Newtonian model 
radius and assuming M = M x , i.e., ip = 1, £ = 0.206666, one obtains £ = 1.12176 and z% — 0.25835. 

B.2.2. Given Radius 

On the other hand , if w e know the true radius, r x = (r, we can determine the mass corresponding to the gravity we 
used by solving Eq. (|B4|) for p and then use Eq. (|B6|) or (|B7|) to determine zi + 1, 

<p = Mi i m = ce (yrrr/ra - 1) . m 

For our parameters we then obtain £ = 1, <p = 0.81440, i.e., M 1 = 1.1401630 M©, z x = 0.22789465. 

B.2.3. Minimal Adjustment 

Alternatively, we could search for a minimum deviation of both p and £ from 1, that is, setting £ = 1/tp in Eq. (|B4[) : 

p 6 + 2C,p 2 -1 = (B9) 

with the physical solution 



p 2 = C/6-4C/C* where C = \J 108 + 4 ^729 + 864C 3 . (B10) 
For our parameters we then obtain £ = l/<p = 1.076353, i.e., z x = 0.246993, M x = 1.30069 M©, and r x = 10.764km. 

B.3. Accretion and Eddington Luminosity 
B.3.1. Eddington Luminosity 

The Eddington luminosity is determined by gravitational acceleration being balanced by radiation pressure on 
electrons: LEdd = ^r 2 gc/ k, where k is the opacity. In Newtonian approximation this computes to 

L Edd =4ncGM/n. (Bll) 

This is also the Eddington luminosity 'at infinity', as there is no redshift in the Newtonian case. In the frame of 
(corrected) stellar surface, the Eddington luminosity taking into account GR gravity is given by 

W,i = (1 + «i) 4kcGM x /k = <p(l + z x ) LEdd • (B12) 

This is the same as the scaling factor £ 2 for any luminosity (Sect. IB.2.1j) . 

B.3.2. Accretion Luminosity 
We summarize the scaling laws for mass, radius, accretion rate and luminosity: 

M x = pM = M x , (B13) 

r x = £r = r 00 /{l + z x ) , r w = (1 + z x ) r x = £ (1 + an) r , (B14) 

L x = ( 2 L = (1 + zif Loo , ioo = L x / (1 + z x f = i 2 Lj (1 + z x f , (B15) 

Mi = £M = (1 + z x ) Mao , = Mi/ (1 + z x ) = i 2 Mj (1 + z x ) . (B16) 



Note that for p> = 1 Eq. (JBTj) leads to M^ = M. 

The accretion luminosity we define by L acc = —M<f>, where M is the accretion rate and <f> the gravitational potential. 
In the Newtonian approximation 4> = —GM/r 2 . For this case we have an accretion luminosity and the ratio of accretion 
luminosity to Eddington luminosity given by 

Lacc = MGM/r , Lacc/Lsdd = Mn/^cr . (B17) 

In GR the gravitational potential is given by 4> = ~c 2 z/{l + z) (from g Iel = d r <j), Eq . iBll lB2l TMisner et alJll973l §25.51. 
Using corrected mass and radius, one obtains 

Lacc,i = M x c 2 z x /(1 + z x ) = Mooc 2 z x = M^c 2 zi , (B18) 

Lacc,l _ MiKCZi _ M oa KCZ 1 _ M£ > 2 KCZ 1 _ Mkcz x (B19) 



L Ed d4 AttGM x (1 + zi) 2 47rGM 00 (l + Zi) iirGMp (1 + z x ) 2 4ttGM (1 + z x ) ' 
where we took advantage of Eq. (|B7[) . For an observer at infinity, the accretion luminosity is reduced by (1 + z x ) 2 : 
Lacc.oo = Lacc,l/ (1 + z x f = M x c 2 z x / (1 + z x f = M^z,/ (1 + z x f = Mpc 2 z x / (1 + zx) 2 , (B20) 
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Thus we obtain the following scaling relations 

-^acc,l -^acc.oo C TZ\ L acc Z\ L &cc (B21) 

i'Edd.l ^Edd.oo GM(1 + Zi) L E dd C(l + z l) L-Edd 

That is, for our example and using ip = 1 the ratio of accretion rate relative to Eddington accretion rate scales by 
0.9934 to the 'GR corrected frame', both at the neutron star surface and in the frame of the observer. 

B.4. Limiting Neutron Star Properties 

Finally, if the entire light curve can be fit accurately enough to observations that both gravity and redshift (z\ ) are 
well determined, we can compute the (gravitational) mass and radius of the neutron star. Using the definitions 

Vi = 1 + z\ and 7 = 1- 1/Vi 2 = 2(ip/£ , (B22) 

we obtain 

<p = y i /(^y/T^) =(V 1 2 -1) 2 /(4V 1 3 C 2 ) =z 1 2 (z 1 + 2) 2 /(4C 2 (z 1 + l) 3 ) , 

£ = 7 / (2Cv^7) = (Vi 2 - 1) /(2CVi) = zi (zi + 2) /(2( ( Zl + 1)) . (B23) 
More generally, mass, M\ and radius, n, are obtained from gravitational acceleration, g, and redshift, z\, by 

M 1 = c 4 z 1 2 (z 1 +2) 2 /(4G 5 (z 1 + l) 3 ) , n =c 2 z 1 (z 1 + 2)/(2 5 (z 1 + l)) . (B24) 

Of course, such a determination would require that all degeneracy of model gravity with accretion rate and metallicity 
would be removed and fitting of the light curve is reliable (and non-degenerate) with respect to nuclear data, opacities, 
equation of state, multidimensional effects, magnetic fields, etc. 



